# This file demonstrates how to calculate the geographic distance
# weights used in Zhukov and Stewart (2013)

# example uses 1998 but this can easily be extended to other settings 
# by adjusting the appropriate parameters


library(maptools)
library(spdep)
library(RANN)
#################################################
## LOAD DATA
#################################################

#pull up our main data file
repdata <- read.table("gw_W.csv", sep=",",header=TRUE)
#subset to a particular year
data <- subset(repdata, repdata$year==1998)
#remove a couple of microstates
data <- data[-which(data$CNTRY_NAME %in% c("Bahrain", "Singapore")),]

## Load map
map <- readShapePoly("shapefiles/d1998")

## Define coordinates
map_crd <- cbind(data$POINT_X, data$POINT_Y) #these are the capitols
map_crd2 <- cbind(data$coor.x, data$coor.y)  #these are the country centroids
IDs <- data$CNTRY_CODE
MIDs <- map$CNTRY_CODE
ifelse(length(unique(MIDs))==length(MIDs),print("NO DOUBLE COUNTS"),print("DOUBLE COUNTS"))



#################################################
## GEOGRAPHIC NEIGHBORS
#################################################


########
## Contiguity
########

## Convert meters to degrees
m2d <- function(m){
out <- (m/1852)/60
return(out)
}
deg <- m2d(500000) ## 500 km

## Create contiguity neighbors (500km snap distance) using the Shape Files
W_nb <- poly2nb(map, row.names =MIDs, queen=TRUE, snap=deg)


#########
## Nearest Neighbors (k = 1)
#########

## Capitals
W_knn1 <- knn2nb(knearneigh(map_crd, k=1), row.names=IDs)
W_knn1_mat <- nb2listw(W_knn1)

## Centroids
W_knn1_2 <- knn2nb(knearneigh(map_crd2, k=1), row.names=IDs)
W_knn1_mat_2 <- nb2listw(W_knn1_2)

#########
## Nearest Neighbors (k = 4)
#########

## Capitals
W_knn4 <- knn2nb(knearneigh(map_crd, k=4), row.names=IDs)
W_knn4_mat <- nb2listw(W_knn4)

## Centroids
W_knn4_2 <- knn2nb(knearneigh(map_crd2, k=4), row.names=IDs)
W_knn4_mat_2 <- nb2listw(W_knn4_2)

#########
## Minimum Distance
#########

## Capitals
dist <- unlist(nbdists(W_knn1, map_crd))
W_dist1 <- dnearneigh(map_crd, d1=0, d2=max(dist), row.names=IDs) 
W_dist1_mat <- nb2listw(W_dist1)

## Centroids
dist_2 <- unlist(nbdists(W_knn1_2, map_crd2))
W_dist1_2 <- dnearneigh(map_crd2, d1=0, d2=max(dist_2), row.names=IDs) 
W_dist1_mat_2 <- nb2listw(W_dist1_2)


#########
## Sphere of influence neighbors
#########

## Capitals
W_del <- tri2nb(map_crd)
W_soi <- graph2nb(soi.graph(W_del, map_crd))
W_soi_mat <- nb2listw(W_soi)



##########################################
#Writing out some Matrices
##########################################


## Contiguity

mat_cont <- nb2mat(W_nb, style="B", zero.policy=T)
mat_cont <- as.matrix(mat_cont)
rownames(mat_cont) <-colnames(mat_cont) <- MIDs
mat_cont

## Min dist

mat_dist <- nb2mat(W_dist1, style="B", zero.policy=T)
mat_dist <- as.matrix(mat_dist)
rownames(mat_dist) <- data$numid
colnames(mat_dist) <- data$numid
mat_dist

## KNN

mat_knn4 <- nb2mat(W_knn4, style="B", zero.policy=T)
mat_knn4 <- as.matrix(mat_knn4)
rownames(mat_knn4) <- data$numid
colnames(mat_knn4) <- data$numid
mat_knn4


## Sphere of Influence
mat_soi <- nb2mat(W_soi, style="B", zero.policy=T)
mat_soi <- as.matrix(mat_soi)
rownames(mat_soi) <- data$numid
colnames(mat_soi) <- data$numid
mat_soi
